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ABSTRACT 

An  entirely  new  approach  to  the  large-eddy  simulation  (LES)  of  high-speed  compressible 
turbulent  flows  is  presented.  In  this  new  approach  to  large-eddy  simulations,  subgrid  scale 
stress  models  go  continuously  to  Reynolds  stress  models  in  the  coarse  mesh/infinite  Reynolds 
number  limit,  where  a  Reynolds  stress  calculation  (RANS)  is  done  in  parallel  with  the  LES  to 
get  an  estimate  of  the  Kolmogorov  length  scale.  Hence,  this  formally  constitutes  a  combined 
LES/time-dependent  RANS  capability.  Here,  the  Reynolds  stress  model  that  is  recovered 
in  the  coarse  mesh  limit  has  a  dependence  on  rotational  strains  through  strain-dependent 
coefficients  as  well  as  through  anisotropic  eddy  viscosity  terms  which  are  usually  neglected 
in  subgrid  scale  models.  Unlike  in  existing  subgrid  scale  models,  this  new  model  has  the 
correct  dependence  on  the  mesh  size  through  the  dimensionless  ratio  of  the  computational 
mesh  size  to  the  Kolmogorov  length  scale.  It  is  this  ratio  that  determines  how  well-resolved 
a  computation  is  in  the  numerical  simulation  of  turbulence  and  should  be  used  to  param¬ 
eterize  subgid  scale  models.  Use  is  made  of  the  Morkovin  hypothesis  in  the  formulation  of 
compressible  subgrid  scale  models.  This  allows  for  the  description  of  supersonic  turbulent 
flows  provided  that  they  are  not  in  the  hypersonic  flow  regime.  A  critical  assessment  of  this 
new  approach  is  provided  along  with  a  discussion  of  traditional  approaches  and  the  prospects 
for  future  applications  to  the  high-speed  compressible  flows  of  technological  importance. 
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1.  INTRODUCTION 


Large-eddy  simulations  of  turbulence  have  been  a  popular  tool  in  turbulence  research 
since  the  ground-breaking  work  of  Smagorinsky  (1963)  three  decades  ago.  In  large-eddy 
simulations  (LES),  the  large,  energy  containing  eddies  are  computed  directly  while  the  small 
scales  are  modeled  since  they  are  believed  to  be  more  universal  in  character  (see  Rogallo 
and  Moin  1984  for  an  interesting  review).  Large-eddy  simulations  were  first  used  in  meteo¬ 
rological  computations  over  thirty  years  ago.  However,  this  was  largely  within  the  context 
of  incompressible  flow.  On  the  other  hand,  the  first  complete  large-eddy  simulation  of  a 
high-speed  compressible  turbulent  flow  was  probably  conducted  by  Erlebacher,  et  al  (1992) 
less  than  a  decade  ago.  This,  however,  was  based  on  a  variable  density  extension  of  the 
Smagorinsky  model  which  was  used  with  the  scale-similarity  model  to  form  a  compressible 
version  of  the  linear  combination  model  (see  Bardina,  Ferziger  and  Reynolds  1983).  Thus, 
these  simulations  were,  to  a  large  extent,  based  on  the  Smagorinsky  model  which  has  a  vari¬ 
ety  of  problems  associated  with  it  (see  Speziale  1997a).  For  example,  the  Smagorinsky  model 
depends  on  the  dimensional  mesh  size  A  so  that  the  subgrid  scale  stress  tensor  r,j  — »  oo 
as  A  oo.  Hence,  a  badly  calibrated  Reynolds  stress  model  -  that  is  overly  dissipative  - 
is  recovered  in  the  coarse  mesh  limit.  Furthermore,  the  Smagorinsky  model  has  no  depen¬ 
dence  on  rotational  strains  and  contains  no  anisotropic  eddy  viscosity  terms.  This  makes 
it  impossible  to  properly  describe  rotating  turbulent  flows;  rotations  substantially  impede 
the  cascade  which  requires  subgrid  scale  models  to  become  substantially  less  dissipative  in 
the  rapid  rotation  limit.  In  addition,  the  anisotropic  eddy  viscosity  terms  are  needed  to 
account  for  backscatter  effects  since  they  are  dispersive  in  character.  Thus,  a  new  approach 
to  LES  is  needed  that  transcends  the  Smagorinsky  model  and  is,  furthermore,  more  suitable 
for  compressible  LES.  This  forms  the  motivation  for  the  present  paper. 

An  entirely  new  approach  to  compressible  large-eddy  simulations  will  be  presented  in  this 
paper.  This  approach  will  require  that  a  Reynolds  stress  calculation  (RANS)  be  conducted 
in  parallel  with  the  LES  in  order  to  get  an  estimate  of  the  Kolmogorov  length  scale  (a  subgrid 
scale  stress  model  will  be  developed  that  depends  on  the  ratio  of  the  computational  mesh  size 
to  the  Kolmogorov  length  scale  which  is  the  parameter  that  determines  how  well  resolved 
a  turbulent  flow  is).  This  only  adds,  at  most,  on  the  average  10%  to  the  computational 
expense.  Since,  the  Reynolds-averaged  turbulent  kinetic  energy  and  dissipation  rate  need 


to  be  obtained  anyhow,  they  are  used  to  partially  parameterize  the  constants  of  the  model 
and  to  non-dimensionalize  the  strain  rates.  Furthermore,  it  should  be  noted  that  since  the 
turbulent  dissipation  rate  is  raised  to  the  1/4  power  in  the  definition  of  the  Kolmogorov  length 
scale,  it  is  quite  feasible  to  get  a  good  estimate  of  this  quantity  with  the  current  generation 
of  Reynolds  stress  models  where  it  only  needs  to  be  calculated  to  within  50%.  The  subgrid 
scale  kinetic  energy  and  dissipation  rate  can  vary  too  much  locally  and,  thus,  are  not  that 
suitable  for  this  purpose.  In  addition,  this  allows  subgrid  scale  stress  models  to  go  more 
continuously  to  Reynolds  stress  models  in  the  coarse  mesh  limit.  The  main  purpose  of  this 
new  approach  is  to  develop  subgrid  scale  stress  models  that  go  continuously  to  state-of-the- 
art  Reynolds  stress  models  in  the  coarse  mesh  limit.  Here,  the  anisotropic  eddy  viscosity  of 
an  explicit  algebraic  stress  model  (see  Gatski  and  Speziale  1993)  will  be  used  for  this  purpose 
(this  formally  constitutes  a  two-equation  model).  The  anisotropic  eddy  viscosity  terms  are 
dispersive  in  character  and  can  account  for  backscatter  effects.  Furthermore,  unlike  in  the 
Smagorinsky  model  or  most  other  subgrid  scale  stress  models,  the  coefficients  of  the  model 
depend  on  both  the  irrotational  and  rotational  strain  rate  invariants  which  allows  for  a  better 
description  of  rotating  flows.  Of  course,  in  the  fine  mesh  limit  as  the  computational  mesh 
size  A  0,  the  subgrid  scale  stress  tensor  r,j  0  -  as  with  existing  subgrid  scale  stress 
models  such  as  the  Smagorinsky  model  -  so  that  a  direct  numerical  simulation  (DNS)  is 
recovered.  However,  this  dependence  is  more  properly  parameterized  by  the  ratio  of  the 
computational  mesh  size  to  the  Kolmogorov  length  scale  as  discussed  above. 

As  far  as  compressibility  effects  are  concerned,  the  Morkovin  hypothesis  will  be  made 
use  of  (see  Morkovin  1964).  In  the  Morkovin  hypothesis,  it  is  assumed  that  compressibility 
effects  only  enter  into  the  description  of  the  turbulence  quantities  through  changes  in  the 
mean  density  (also  see  Cebeci  and  Smith  1974).  The  Morkovin  hypothesis  has  been  validated 
in  wall-bounded  turbulent  flows,  provided  they  are  not  in  the  hypersonic  flow  regime  where 
the  external  Mach  number  Maoo  >  5.  For  the  flat  plate  turbulent  boundary  layer,  it  is 
approximately  valid  for  external  Mach  numbers  Mdoo  <  8  (see  Zhang,  So,  Speziale  and 
Lai  1993).  Since  the  external  Mach  number  is  somewhat  less  at  the  subgrid  scale  level, 
the  Morkovin  hypothesis  is  probably  a  good  approximation  for  the  development  of  subgrid 
scale  models,  in  general  turbulent  flows,  for  external  Mach  numbers  Ma^o  <  8.  Thus,  the 
supersonic  flow  regime  -  as  well  as  the  beginning  part  of  the  hypersonic  flow  regime  -  can  be 


well  described  by  this  hypothesis.  Of  course,  this  hypothesis  breaks  down  entirely  in  strongly 
hypersonic  flows  where  real  gas  effects  and  dissociation  must  be  accounted  for.  Turbulent 
dilatational  effects  -  which  are  neglected  in  the  Morkovin  hypothesis  -  must  also  be  accounted 
for  in  compressible  turbulent  flows  where  the  external  Mach  number  Maoo  is  greater  than 
4  or  5  (more  precisely,  when  the  turbulence  Mach  number,  which  will  be  discussed  later,  is 
greater  than  approximately  0.3).  While  the  inclusion  of  the  dilatational  dissipation  can  be 
useful  in  free  turbulent  shear  flows  for  large  Mach  numbers  (see  Sarkar,  Erlebacher,  Hussaini 
and  Kreiss  1991  and  Zeman  1990),  it  is  counterproductive  in  wall-bounded  turbulent  flows 
where  it  can  cause  a  degradation  of  the  predictions.  Thus,  turbulent  dilatational  models  will 
not  be  made  fuU  use  of  here  for  the  development  of  subgrid  scale  models. 

The  subgrid  scale  modeling  for  this  combined  LES /time-dependent  RANS  capability  will 
be  discussed  in  detail  in  the  sections  to  follow.  It  has  the  potential  to  bridge  the  gap  between 
DNS,  LES  and  RANS.  A  fuU  discussion  of  the  implications  for  turbulence  research  will  be 
provided  in  the  last  section. 


2.  A  REVIEW  OF  COMPRESSIBLE  LARGE-EDDY  SIMULATIONS 


The  full  equations  of  motion  for  an  ideal  gas  will  be  considered  (cf.  Cebeci  and  Smith 

1974): 


Continuity 


^  =  0 


(1) 


Momentum 


Q 

+  {pUiUj),j  =  -P,i  +  (TijJ 

2 

K-  =  +  Mi,,)] 


(2) 


Energy 


dt 


{pCvT)  +  (puiCvT)^i  =  —pui^i  +  $  + 


[$  =  O-ijUij  =  --p{ui,i)^  +  p{Uij  +  Uj,i)uij] 


(3) 


State 


p  =  pRT 


(4) 


where 


p  =  mass  density 
Ui  =  velocity  vector 
p  =  thermodynamic  pressure 
p,  =  dynamic  viscosity 
(Tij  =  viscous  stress  tensor 
T  =  absolute  temperature 
K  =  thermal  conductivity 
R  =  ideal  gas  constant 
Cv  =  specific  heat  at  constant  volume 
$  =  viscous  dissipation  function 

(  ),.•  ^  *(  )• 

For  any  fiow  variable  J-,  we  can  introduce  the  decomposition 


T  =  7  +  r 


(5) 


5 


where  ^  is  a  standard  filter  and  a  prime  represents  a  subgrid  scale  fluctuation.  These, 
respectively,  represent  the  large  and  small  scale  fields  of  the  turbulence.  Alternatively,  the 
decomposition 

T  =  (6) 

can  be  introduced  where  _ 

(7) 

P 

is  the  Favre  (or  mass-weighted)  filter.  Here, 

^  ^  0,  ^  ^  0 


and 

f  54  0,  0, 

in  general.  A  filtered  quantity  is  given  by 

T  =  I  (3(x-x*,A):F(x*)dV  (8) 

Jd 

In  (8),  A  is  the  computational  mesh  size  and  G  is  a  filter  function  which  is  normalized  as 
follows: 

f  G(x-x*,A)dV  =  1.  (9) 

J  D 

This  guarantees  that  G  becomes  a  Dirac  delta  sequence  in  the  limit  as  A  — >  0: 


Hm  G(x  -  X*,  A)«^(x*)d"x*  = 


where  5(x  —  x*)  is  the  Dirac  delta  function.  Direct  numerical  simulations  (DNS)  are,  thus, 
recovered  in  the  fine  mesh  limit.  Due  to  the  Riemann-Lebesgue  Theorem,  (8)  substantially 
reduces  the  amplitude  of  the  high-wavenumber  Fourier  components  in  space  of  any  flow 
variable  (consequently,  J-  represents  the  large  scale  part  of  The  filter  function  G 
has  usually  been  taken  to  be  a  Gaussian  filter  in  infinite  domains  or  a  piecewise  continuous 
distribution  of  bounded  support  in  compact  domains  (in  the  latter  case,  the  simple  box  filter 
has  been  commonly  used  with  finite  difference  methods;  see  Deardorff  1970).  The  box  filter 
on  a  non-uniform  mesh  is  given  by 


G(x  -  X*, 


1/8  A^ 

0, 


I®.-  -  ®i*l  <  Axi 

i®i  -  ®I*1  >  Ax; 


(10) 


for  t  =  1,2,3  where  A*,.,  is  A^;,  Aj,  and  A^,  respectively,  and  where  A  is  given  by 

A  =  (A.A„A,)''" 

(Ax,  Ay  and  Ax  are  the  mesh  sizes  in  the  x,  y  and  z  directions,  respectively,  obtained  after 
a  coordinate  transformation). 

The  filtered  continuity,  momentum  and  energy  equations  take  the  form 


Continuity 


Momentum 


Energy 


%  +  “ 


Q 

+  {pUiUj),j  =  -p,,-  +  O'ijj  -  {pTij),j 

(r,j  =  u'lu'f,  p  =  pRT) 


Q  ^  _ 

— (pC„r)  +  {puiC„f),i  =  -pui,i  -  pu’ii 


(11) 


(12) 


(13) 


+  ^  +  ('^^,«),«  Qi,i 
(g.  =  pau?T") 

where  turbulent  fluctuations  in  Cv  have  been  neglected.  Here,  Tij  is  the  Favre-filtered  subgrid 
scale  stress  tensor  and  Qi  is  the  Favre-filtered  subgrid  scale  Reynolds  heat  flux;  $  is  the 
filtered  dissipation  function.  The  filtered  dissipation  function  is  given  by: 


$  =  <TijUi,j  +  ar'iju'i 


(14) 


(TijUiJ  +  +  P£ 


where  e  =  (r'lju'i  jj p  is  the  filtered  turbulent  dissipation  rate.  Hence,  $  — >  pe  as  the  Reynolds 
number  Re  oo. 

In  high-Reynolds-number  turbulent  flows,  the  molecular  diffusion  terms  are  dominated 
by  the  turbulent  transport  terms  except  in  a  thin  sublayer  near  walls.  If  we  assume  in  this 
region  that  fluctuations  in  the  viscosity,  thermal  conductivity  and  density  can  be  neglected 
we  can  then  make  the  approximations: 


o-ti  =  —^puk,kSij  +  fi{uij+Uj,i) 
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Thus,  in  order  to  achieve  closure,  models  are  needed  for: 


(1)  The  Favre-filtered  subgrid  scale  Reynolds  stress,  T,j 

(2)  The  Favre-filtered  subgrid  scale  Reynolds  heat  flux,  Qi 

(3)  The  subgrid  scale  mass  flux,  u" 


(4)  The  subgrid  scale  dissipation  rate,  e 

(5)  The  subgid  scale  pressure-dilatation  correlation, 

This  makes  the  problem  of  conducting  compressible  large-eddy  simulations  far  more  difficult 
than  its  incompressible  counterpart. 

The  full  form  of  the  subgrid  scale  stress  tensor  is  given  by 

Tij  —  Lij  -|-  Cij  -f  Rij  (15) 

where 

Lij  =  UiUj  -  UiUj  (16) 

Cij  =  tijtij  +  u'/uj  (17) 

Rij  s  (18) 

are,  respectively,  the  Leonard  stresses,  subgrid  scale  cross  stresses  and  subgrid  scale  Reynolds 
stresses  (see  Leonard  1974  and  Ferziger  1976). 

The  Smagorinsky  model  has  been  used  for  either  the  entire  deviatoric  subgrid  scale  stress 
tensor  or  for  the  deviatoric  part  of  the  subgrid  scale  Reynolds  stress  tensor.  In  the  latter 
case,  for  compressible  flows,  the  deviatoric  part  of  the  Reynolds  subgrid  scale  stress  tensor 
is  modeled  as 

oRij  =  -2C^,A\2SklSMY^\SiJ  -  ^SmmSij)  (19) 

where 

1  /  du{  dui\ 


is  the  rate  of  strain  tensor.  In  the  incompressible  limit,  Sij  — >  Sij  and,  furthermore. 


Smm  =  V  •  u  =  0 


(20) 


SO  that  the  incompressible  form  of  the  Smagorinsky  model  is  recovered.  Here,  (7,  is  the 
Smagorinsky  constant  that  assumes  the  approximate  value  of  0.1  in  many  flows  but  can  vary 
by  as  much  as  a  factor  of  three.  When  the  Smagorinsky  model  is  used  for  the  subgrid  scale 
Reynolds  stress  tensor  it  is  usually  used  in  conjunction  with  the  scale  similarity  model  for 
the  subgrid  scale  Leonard  and  cross  stresses  given  by  (see  Bardina,  Ferziger  and  Reynolds 
1983) 

Lij  +  Cij  =  UiUj  -  UiUj  (21) 


whose  coefficient  is  unity  in  order  to  preserve  Galilean  invariance  (see  Speziale  1985).  Birin- 
gen  and  Reynolds  (1981)  and  Moin  and  Kim  (1982)  had  violated  Galilean  invariance  by 
applying  the  Smagorinsky  model  to  the  sum  of  the  cross  stresses  and  Reynolds  subgrid  scale 
stresses  (see  Speziale  1985). 

The  subgrid  scale  Reynolds  stress  tensor  can  be  split  into  isotropic  and  deviatoric  parts 
as  follows: 


Rij  —  ^Rmm^ij  "b  oRij 


where  Sij  is  again  the  Kronecker  delta.  The  isotropic  part,  Rmm  =  the  Reynolds 

subgrid  scale  stress  tensor  can  be  absorbed  into  the  pressure.  In  many  flows,  it  is  an  order 
of  magnitude  smaller  than  the  thermodynamic  pressure  so  it  has  been  neglected  in  many 
applications.  During  the  last  decade,  Yohizawa  (1986)  developed  a  model  for  Rmm  which  is 
given  by 

Rmm  =  -2Clp^^SmnSmn  (22) 


where  Cj  is  a  constant.  Speziale,  Erlebacher,  Zang  and  Hussaini  (1988)  showed  that  this 
proposed  model  performs  poorly.  A  scatterplot  of  this  model  compared  with  DNS  is  shown 
in  Figure  1  showing  the  poor  level  of  correlation. 

The  Smagorinsky  model  has  several  deficiencies  that  can  be  summarized  as  follows: 

(1)  The  Smagorinsky  constant  is  not  in  reality  a  constant.  It  can  vary  by  as  much 
as  a  factor  of  two  or  three  feom  flow  to  flow.  This  is  because  the  Smagorinsky  model 
is  badly  parameterized  (see  Speziale  1997a).  Furthermore,  it  only  correlates  with  DNS 


at  the  50%  level.  To  get  an  idea  of  how  poor  this  result  is,  the  correlation  between  the 
functions  y  —  x  and  y  =  c  ^  on  the  interval  [0,  1]  is  more  than  50%  despite  the  fact  that 
they  are  qualitatively  different  functions  (one  is  monotonically  increasing  while  the  other  is 
monotonically  decreasing)! 


(2)  The  Smagorinsky  model  does  not  depend  on  rotational  strains  through  the  invariant 
^  a  (Wij  =  Kduifdxj-dujidxi)  is  the  Favre-filtered  vorticity  tensor)  and,  fur¬ 

thermore,  has  the  wrong  dependence  on  the  irrotational  strain  rate  invariant  7^  oc 
For  Reynolds  stress  models  in  equilibrium,  the  eddy  viscosity  reduces  to 


vj  oc 


3 

3-27/2-1-  6^2 


(see  Gatski  and  Speziale  1993). 

(3)  The  dependence  on  the  computational  mesh  size  A  should  be  through  the  dimen¬ 
sionless  ratio  AfLK-  What  determines  how  well  a  computation  is  resolved  in  the  numerical 
simulation  of  turbulence  is  whether  or  not  the  grid  size  is  small  ( or  large)  compared  to  the 
Kolmogorov  length  scale.  The  dimensional  dependence  on  A  in  the  Smagorinsky  model  is 
simply  incorrect.  In  the  Smagorinsky  model,  t,j  — >•  oo  as  A  — >  oo.  Hence,  a  badly  calibrated 
Reynolds  stress  model  is  obtained  in  the  coarse  mesh  limit.  The  model  becomes  far  too 
dissipative  as  the  mesh  becomes  coarse. 

In  so  far  as  point  (2)  is  concerned,  this  makes  it  impossible  for  the  Smagorinsky  model 
to  properly  describe  rotating  flows.  For  example,  it  is  well  known  that  in  an  incompressible 
rapidly  rotating  isotropic  turbulence,  the  cascade  is  essentially  shut  off  so  that  the  turbulence 
undergoes  a  linearly  viscous  decay  (see  Speziale,  Mansour  and  Rogallo  1987).  Hence,  it  is 
possible  to  conduct  direct  simulations  even  at  high  turbulence  Reynolds  numbers.  The 
Smagorinsky  model  is  far  too  dissipative  in  this  case  where  it  can  yield  results  that  are 
completely  erroneous.  For  a  rapidly  rotating  isotropic  turbulence,  the  Smagorinsky  constant 
is  essentially  zero  except,  perhaps,  at  astronomically  high  Reynolds  numbers  or  for  extremely 
coarse  meshes. 

In  the  next  section,  a  new  methodology  for  the  large-eddy  simulation  of  compressible 
turbulent  flows  will  be  presented.  There  is  no  question  that  a  new  approach  that  transcends 
the  Smagorinsky  model  is  needed.  Even  the  dynamic  subgrid  scale  model  of  Germano, 
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Piomelli,  Moin  and  Cabot  (1991)  does  not  overcome  these  shortcomings  as  discussed  in 
Speziale  (1997a). 
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3.  A  NEW  APPROACH  TO  COMPRESSIBLE  LARGE-EDDY  SIMULATIONS 


The  new  approach  to  compressible  LES  that  is  being  proposed  has  subgrid  scale  stress 
models  that  are  of  the  following  form: 


Tij  =  [1  -  exp{-/3A/LK)T 


^KSy  -  01/(1;, 0^(4-  - 


+  Ti 


tj 


(23) 


where  r-j  represents  the  anisotropic  part  of  the  subgrid  scale  stress  tensor.  Here,  an  overtilde 
represents  a  Favre-filter  whereas 

,  =  o,(44)V^|,  (  =  (24) 

where  Sij  and  Wij  are  the  Favre-filtered  rate  of  strain  and  vorticity  tensors,  A  is  the  com¬ 
putational  mesh  size,  Lk  =  (where  V  =  Jtjp)  is  the  Kolmogorov  length  scale, 

and  13,  n,  aj,  ©2  aiid  0:3  are  constants  (ai,  02  and  03  are  obtained  from  a  Reynolds  stress 
model  along  with  the  function  /).  Here,  K  and  £  represent  the  Reynolds-averaged  turbu¬ 
lent  kinetic  energy  and  dissipation  rate  obtained  from  a  Reynolds  stress  calculation  with  the 
two-equation  model  equivalent  to  that  given  above  in  the  coarse  mesh  limit  as  AJLk  —*■  00. 

In  the  coarse  mesh  limit,  a  Reynolds  stress  model  given  by 

Ilf  =  IkSh  -  •*  (26) 

is  recovered  which  is  an  explicit  algebraic  stress  model  (see  Gatski  and  Speziale  1993  and 
Speziale  1996  for  the  compressible  extension  of  this  model).  The  turbulent  dissipation  rate 
£  -  and,  hence,  the  turbulent  kinetic  energy  K  -  have  to  be  obtained  anyway  in  order  to 
get  an  estimate  of  the  Kolmogorov  length  scale  Lk.  Since,  the  Kolmogorov  length  scale 
Lk  =  the  dissipation  rate  only  has  to  be  estimated  to  within  50%  with  the 

modeled  dissipation  rate  equation  to  get  a  good  estimate  of  the  Kolmogorov  length  scale 
(the  dissipation  rate  is  raised  to  the  1/4  power  as  mentioned  before).  This  is  quite  feasible 
with  state-of-the-art  Reynolds  stress  models.  Thus,  this  methodology  requires  that  a  RANS 
calculation  be  done  in  parallel  with  the  LES.  This  will,  in  most  circumstances,  only  add 
at  most  10%  to  the  computational  expense.  Here,  we  parameterize  the  model  in  terms  of 
the  Reynolds- averaged  turbulent  kinetic  energy  and  dissipation  rate  since  the  subgrid  scale 
turbulent  kinetic  energy  and  dissipation  rate  can  vary  too  much  locally.  This  model  has 


been  written  before  in  the  shorthand  notation  as  (see  Speziale  1997b) 

Tij  =  [1  -  ea!p(-/3A/iA-)]"4^^  (26) 

where  is  a  Reynolds  stress  model  that  is  written  partially  in  terms  of  filtered  fields.  An 
explicit  algebraic  stress  model  is  used  for  this  purpose  as  discussed  above. 

The  anisotropic  eddy  viscosity  terms  take  the  form 

=  [1  -  +  Wit'Sti) 

1 

where  again  the  overtilde  represents  a  Favre-filtered  quantity  whereas  K  and  6  are  the 
Reynolds- averaged  turbulent  kinetic  energy  and  dissipation  rate  obtained  from  a  Reynolds 
stress  calculation  (a4  and  as  are  constants).  In  the  coarse  mesh  limit  as  A./Lk  -»  oo,  the 
anisotropic  eddy  viscosity  terms  of  an  explicit  algebraic  stress  model  are  recovered,  extended 
to  compressible  flows  (see  Gatski  and  Speziale  1993  and  Speziale  1996).  These  terms  are 
dispersive  in  character  and  can  account  for  backscatter  effects  (cf.,  Clark,  Ferziger  and 
Reynolds  1979) 

For  incompressible  Reynolds  stress  models  in  equilibrium 

A  singularity  can  occur  when  this  expression  is  applied  to  turbulent  flows  where  there  are 
significant  departures  from  equilibrium.  Gatski  and  Speziale  (1993)  introduced  the  simple 
regularization 

- ? - _ - -  (29) 

3-27/2-1-6^2  3  +  ^2  +  6^2^2  +  6^2 

which  is  obtained  by  a  Taylor  series  expansion.  For  turbulent  flows  in  equilibrium  where 
7/,^  <  1,  it  yields  results  that  are  indistinguishable  from  (28)  where  it  formally  applies. 
But  it  is  regular  and  computable  for  all  values  of  7/  and  More  recently,  Speziale  and  Xu 
(1996)  obtained  expressions  via  a  formal  Fade*  approximation  that  builds  in  some  limited 
agreement  with  the  Rapid  Distortion  Theory  (RDT)  solutions  for  plane  shear  and  plane 
strain  turbulence.  The  constants  in  this  Reynolds  stress  model  are  given  by 

ai  =  0.374,  a2  =  0.145,  aj  =  0.308,  (30) 
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04  =  0.115,  05  =  0.108.  (31) 

These  constants  are  based  on  benchmark  flows  in  the  incompressible  limit  since  the  Morkovin 
hypothesis  is  being  applied  where  compressible  models  are  a  variable  density  extension  of 
their  incompressible  counterparts. 

The  Reynolds- averaged  turbulent  kinetic  energy  K  and  dissipation  rate  E  are  obtained 
from  modeled  versions  of  their  transport  equations  which  take  the  compressible  form  (cf. 
Speziale  1996  and  Speziale  and  Sarkar  1991): 

+[{p+ 

-i?£i,„  -0£#|:+  [()r+^)£,,]  ,  (33) 

where  V  =  -T^pduijdxj  is  the  turbulence  production,  =  C^pK^IE  is  the  eddy  viscosity 
and  Csi,  Cs2,  and  ae  are  constants  that  assume  the  values  of  0.09,  1.44,  1.83,  1.0 
and  1.3,  respectively.  These  equations  have  served  as  a  cornerstone  for  two-equation  models 
in  the  incompressible  limit.  In  order  to  integrate  this  model  to  a  wall  it  is  oiJy  necessary  to 
remove  the  singularity  in  the  destruction  term  that  appears  on  the  right-hand-side  of  (33) 
with  the  coefficient  C£2  (see  Speziale  and  Abid  1995).  No  ad  hoc  wall  damping  functions  are 
needed  in  the  Reynolds  stress  model  This  is  accomplished  by  replacing  Cs2  with  the  term 

Cc2[l  -  e®p(-Ri,/10)]  (34) 

where  Ry  =  K^/'^ylv  given  that  y  is  the  coordinate  normal  to  the  wall.  In  many  applications, 
a  small  vortex  stretching  term  has  been  added  to  (33)  to  make  the  calculations  better  be¬ 
haved.  It  removes  the  singularity  in  plane  stagnation  point  turbulent  flows  and,  furthermore, 
allows  for  the  description  of  both  the  log-layer  and  homogeneous  turbulence  in  equilibrium 
with  a  simple  unifled  model  where  it  is  not  necessary  to  solve  the  cubic  equation  arising  out 
of  the  consisistency  condition  (see  Abid  and  Speziale  1996  and  Speziale,  Jongen  and  Gatski 
1997).  The  dilatational  dissipation  has  been  neglected  since  it  leads  to  a  degradation  of  the 
results  in  wall-bounded  flows  (see  Speziale  1996).  The  dilatational  dissipation  was  proposed 
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by  Zeman  (1990)  and  Sarkar,  Erlebacher,  Hussaini  and  Kreiss  (1991)  largely  based  on  free 
turbulent  shear  flows. 

It  is  worthwhile  to  note  that  these  equations  are  consistent  with  the  limit  of  compressed 
isotropic  turbulence.  For  this  problem,  the  mean  velocity  gradient  tensor  is  given  by 


/  |F  0  0 

'“«.i  ~  1  ®  3^  ® 

V  0  0  |F 

\ 

/ 

(35) 

where  F  is  the  expansion/compression  rate  which  is  constant.  The  model  provided  herein 

reduces  approximately  to  the  simple  coupled  ODE’s: 

2 

if  =  -jrjr 

(36) 

II 

1 

CO) 

(37) 

for  \£\KolSo  ^  1.  The  short-time  solution  to  Eqs.  (36) 

-  (37)  is  given  by: 

K  =  Koexp 

\  0  ) 

1 

(38) 

E  =  5oexp 

(39) 

A  =  Aoexp 

(40) 

where  A  =  jS  is  the  integral  length  scale.  These  are  identical  to  the  results  obtained 
by  Reynolds  (1987)  based  on  Rapid  Distortion  Theory  (RDT).  In  contrast  to  these  results,  a 
variable  density  extension  of  the  commonly  used  second-order  closures  erroneously  predicts 
that 

A  =  Aoexp(— 0.04rt)  (41) 

(see  Reynolds  1987).  According  to  (41),  the  integral  length  scale  will  decrease  under  an 
expansion  (F  >  0)  and  increase  under  a  compression  (F  <  0)  —  results  that  are  clearly  in 
error  as  first  pointed  out  by  Reynolds  (1987). 

It  is  also  worth  noting  that  in  a  rapidly  rotating  incompressible  isotropic  turbulence, 
f{7},0  — >  0  so  Tij  0  yielding  a  DNS.  This  results  from  the  fact  that 


in  rotating  frames  where  is  the  rotation  rate  of  the  frame  and  emji  is  the  permutation 
tensor  (hence,  ^  ~  in  a  rapidly  rotating  flow  with  angular  velocity  il).  As  mentioned 
earlier,  in  a  rapidly  rotating  incompressible  isotropic  turbulence  the  energy  cascade  is  es¬ 
sentially  shut  off  so  that  direct  numerical  simulations  can  be  conducted  with  a  128^  mesh 
even  at  high  turbulence  Reynolds  numbers.  A  128^  mesh  forms  a  cornerstone  of  this  new 
approach  to  LES  as  will  soon  be  discussed.  In  contrast,  the  Smagorinsky  model  is  far  too 
dissipative  so  it  yields  incorrect  results  for  this  problem. 

The  grid  function 

[1  -  exp{-pA/LK)r  (43) 

bridges  the  gap  between  DNS,  LES  and  RANS  where  Lk  is  the  Kolmogorov  length  scale 
Lk  =  estimated  from  a  Reynolds  stress  calculation  (again,  A  is  the  computational 

mesh  size).  In  the  limit  as  A/Lk  oo  the  grid  function  goes  to  one  and  we  recover  a 
Reynolds  stress  model  whereas  in  the  limit  as  A/Lk  0,  it  goes  to  zero  and  we  formally 
recover  a  DNS.  Actually,  when  AfLx  is  of  order  one,  we  should  have  a  DNS  (this  has  been 
built  into  the  calibration).  Since  Lk  =  jS  where  Ri  =  juS  is  the  turbulence 

Reynolds  number,  AJLk  — >  oo  as  Rf  — >  oo  (thus,  we  recover  a  Reynolds  stress  model  in  the 
coarse  mesh/infinite- Reynolds-number  limit).  For  the  initial  calculations,  n  has  been  taken 
to  be  one  and  ^  has  been  calibrated  as  follows: 

P  «  0.001.  (44) 

A  power  law  for  the  grid  function  has  been  theoretically  obtained  using  Renormalization 
Group  methods  (Woodruff  and  Hussaini,  Private  Communication).  Note  that  for  AJLk  < 
100,  we  approximately  obtain  a  power  law  from  (43)  via  a  Taylor  expansion.  Most  practical 
LES  are  conducted  for  A/Lk  =  10  —  100. 

The  modeling  of  the  other  higher-order  correlations  that  are  needed  for  the  closure  of  the 
compressible  equations  is  along  the  same  lines  as  that  for  the  Reynolds  stress  models.  In 
effect,  for  this  approach,  a  subgrid  scale  model  is  nothing  more  than  a  scaled  down  version 
of  the  Reynolds  stress  model  -  scaled  down  by  the  grid  function  to  have  less  dissipation.  For 
the  subgrid  scale  heat  flux 

f,-  (45) 

I'tt  e 


16 


The  turbulent  Prandtl  number  Prj  ^  0.9.  In  the  coarse  mesh  limit,  as  the  grid  function 
goes  to  one,  the  standard  gradient  transport  model  for  the  Reynolds  heat  flux  is  recovered. 
The  same  is  true  of  the  subgid  scale  mass  flux  which  is  modeled  as 

<  =  [1  -  («) 

(see  Cebeci  and  Smith  1974).  Here,  the  constant  CTp  is  approximately  equal  to  0.5.  The 
subgrid  scale  dissipation  e  and  Reynolds-averaged  dissipation  E  are  related  by  the  grid 
function  as  follows 

£  =  [1  _  exp{-fiA/LK)rE.  (47) 

The  pressure-dilatation  correlation  is  modeled  as 

=  [1  —  exp{—(3A/LK)y'{<i2pTijU{jMt  +  aspeMf}  (48) 

where 

a2  =  0.15,  aa  =  0.2  (49) 

and 

Mt  =  {2Kl'fRfy^^  (50) 

is  the  turbulence  Mach  number  (7  is  the  ratio  of  speciflc  heats  and  R  is  the  ideal  gas 
constant).  In  the  coarse  mesh  limit,  as  the  grid  function  goes  to  one,  a  recently  proposed 
Reynolds-averaged  model  is  obtained  for  the  pressure-dilatation  correlation  that  is  suitable 
for  compressible  shear  flows  (see  Sarkar  1992).  Thus,  we  have  a  complete  closure  to  the 
compressible  filtered  equations  of  motion.  The  full  Reynolds  stress  model  that  is  recovered 
in  the  coarse  mesh  limit  has  been  tested  successfully  in  a  variety  of  benchmark  compressible 
flows  (see  Speziale  1996). 

Some  comments  are  needed  concerning  the  choice  of  a  filter  in  this  new  approach  to 
large-eddy  simulations  and  the  melding  together  of  spatial  filtering  in  LES  and  Reynolds 
averaging  in  RANS.  We  want  a  filter  that  yields  the  minimum  contamination  of  the  large 
scales.  The  reason  for  this  is  clear;  defiltering  must  be  avoided  since  it  constitutes  an  ill- 
posed  mathematical  problem  (see  Speziale  1997b).  The  purpose  of  practical  LES  is  to  predict 
the  Reynolds- averaged  fields.  In  order  to  do  so,  the  filtered  velocity,  which  is  calculated, 
must  invariably  be  used  to  estimate  the  large-scale  part  of  the  instantaneous  velocity  which 


then  yields  the  Reynolds-averaged  fields  through  appropriate  ensemble  or  time  averages. 
The  filtered  equations  of  motion  (11)-(13)  are  of  the  same  form  as  the  Reynolds-averaged 
equations.  In  the  coarse  mesh  limit,  the  ramp  function  will  be  one  and  the  model  will  be  so 
dissipative  that  a  RANS  calculation  will  be  automatically  recovered  with  a  state-of-the-art 
Reynolds  stress  model.  It  is  envisioned  that  ensemble  averages  will  be  taken  even  if  we  are 
conducting  a  time- dependent  RANS.  Thus,  we  do  not  need  to  know  the  effect  of  the  filter  - 
which  can  never  be  fuUy  known  in  complex  geometries  -  except,  perhaps  for  model  calibration 
in  benchmark  flows.  This  allows  us  to  meld  together  the  LES  and  RANS  methodologies  which 
are  normally  treated  as  disparate  approaches.  In  both  of  these  approaches  we  calculate  what 
is  tantamount  to  the  large-scale  velocity  field  -  through  the  same  basic  equations  of  motion 
-  and  then  obtain  the  Reynolds-averaged  fields  through  ensemble  averages  (time  averages  in 
a  statistically  steady  turbulence).  The  large  scales  make  the  dominant  contribution  to  the 
most  pertinent  fields  such  as  the  turbulent  kinetic  energy.  A  minimum  contamination  of  the 
large  scales  can  be  accomplished  with,  of  the  order  of,  a  128^  computational  mesh  using  a 
filter  with  a  compact  support  -  such  as  the  box  filter  -  which  has  a  small  filter  width  of,  for 
example,  two  mesh  points.  Some  of  the  previously  conducted  coarse  grid  LES  (which  has 
typically  had  no  more  than  32^  mesh  points)  must  be  avoided  wherein  the  filter  width  has, 
at  times,  been  as  much  as  25%  of  the  computational  domain,  significantly  contaminating  the 
large  scales.  Besides,  recent  increases  in  computational  capacity  have  begun  to  make  128^ 
computations  much  more  feasible  for  engineering  calculations  (a  small  compromise  to  100^ 
computations  can  always  be  made).  In  addition,  it  should  be  noted  that  practical  LES  -  in 
complex  geometries  -  will  require  the  use  of  finite  difference  techniques  with  a  compact  filter 
where  we  will  never  make  explicit  use  of  the  filter.  These  finite  difference  methods  should, 
furthermore,  be  based  on  fourth-order  accurate  finite  difference  schemes  for  better  accuracy. 
Spectral  methods  have  to  be  abandoned  if  complex  turbulent  flows  are  to  be  addressed.  With 
this  new  methodology,  the  gap  between  DNS,  LES  and  RANS  can  be  bridged  (see  Figure 


4.  CONCLUSION 


A  new  approach  to  the  large-eddy  simulation  of  the  high-speed  compressible  turbulent 
flows  of  technological  importance  has  been  proposed.  In  this  new  methodology,  subgrid 
scale  models  go  continuously  to  Reynolds  stress  models  in  the  coarse  mesh/infinite  Reynolds 
number  limit.  Hence,  with  this  new  methodology  it  is  possible  to  achieve  the  long  held  dream 
of  going  continuously  from  a  large-eddy  simulation  to  a  Reynolds  stress  calculation  as  the 
mesh  becomes  coarse  or  the  Reynolds  number  becomes  extremely  large.  It  is  firmly  believed 
that  in  complex  wall-bounded  turbulent  flows  -  especially  with  flow  separation  where  wall 
functions  cannot  be  used  —  the  best  that  one  can  do,  at  this  time,  for  the  extremely  high 
Reynolds  numbers  encountered  in  many  technological  applications  (such  as  some  naval  flows 
where  Re  ~  0(10^))  is  a  RANS  computation  since  the  crucial  wall  layer  cannot  be  resolved. 
At  more  moderate  Reynolds  numbers,  LES  is  possible  (of  course,  at  low  turbulence  Reynolds 
numbers  it  is  possible  to  conduct  a  DNS).  The  RANS  model  that  was  presented  in  this  study 
is  a  two-equation  model  that  contains  some  of  the  most  recent  developments  in  compressible 
turbulence  modeling.  In  the  incompressible  limit  it  coUapses  to  an  explicit  algebraic  stress 
model  which  is  a  two-equation  model  that  is  consistent  with  a  state-of-the-art  second-order 
closure  model  in  the  limit  of  homogeneous  turbulence  in  equilibrium.  A  new  approach  is 
unquestionably  needed  that  transcends  the  Smagorinsky  model.  The  Smagorinsky  model  has 
a  variety  of  deficiencies  as  outlined  in  this  paper  (also  see  Speziale  1997a).  The  only  reason 
to  believe  that  it  has  been  successful  in  previous  applications  is  because  it  drains  enough 
energy  to  approximately  account  for  the  energy  cascade  to  the  scales  that  are  left  unresolved 
-  an  effect  that  is  achieved  by  the  ad  hoc  adjustment  of  the  Smagorinsky  constant. 

Some  brief  remarks  are  in  order  concerning  the  role  of  direct  and  large-eddy  simulations 
in  turbulence.  There  is  no  question  that  DNS  —  and  the  computer  in  general  —  has  revolu¬ 
tionized  the  study  of  turbulence.  DNS  has  already  shed  new  light  on  the  physics  of  a  range 
of  basic  turbulent  flows  and  the  future  potential  is  enormous.  It  already  appears  that  in 
the  not  too  distant  future,  DNS  will  entirely  replace  basic  benchmark  physical  experiments 
for  homogeneous  turbulence,  near-wall  turbulent  flows  and  basic  compressible  flows  at  lower 
turbulence  Reynolds  numbers.  However,  it  appears  that  DNS  will,  for  a  long  time  to  come, 
be  limited  to  relatively  simple  geometries  and  low  to  moderate  turbulence  Reynolds  num¬ 
bers.  Direct  simulations  of  the  complex  turbulent  flows  of  technological  importance,  at  high 


turbulence  Reynolds  numbers,  could  require  the  generation  of  data  bases  with  upwards  of 
10^°  numbers.  Thus,  it  is  crucial  that  large-eddy  simulations  be  made  to  work.  As  far  as 
LES  is  concerned,  it  must  be  said  that  it  has  never  lived  up  to  its  initial  promise.  The  way 
that  traditional  LES  has  been  formulated  is  not  predictive  and  is  probably  only  suitable  for 
doing  less  expensive  parametric  studies  of  benchmark  direct  simulations  once  the  reliability 
of  the  subgrid  scale  model  has  been  established  by  DNS  for  the  baseline  case.  In  order  to 
solve  the  complex  turbulent  flows  of  technological  importance,  an  entirely  new  approach  to 
LES  is  direly  needed. 

Finally,  some  comments  are  warranted  concerning  supersonic  turbulent  flows  -  partic¬ 
ularly  in  wall-bounded  geometries.  Since  there  is  considerable  heating  in  these  flows,  the 
effective  Reynolds  number  (through  a  rise  in  the  kinematic  viscosity)  is  not  extremely  large 
even  in  many  of  the  high-speed  compressible  flows  of  technological  importance.  Thus,  the 
prospects  for  being  able  to  conduct  a  large-eddy  simulation  in  many  practical  flow  situations 
are  fairly  good  for  compressible  turbulence.  This  would  include  the  flow  around  high-speed 
aircraft  and  missiles  to  name  just  a  few  applications.  For  those  cases  where  LES  is  not 
feasible,  a  time-dependent  RANS  can  be  always  be  conducted  which  is  expected  to  be  far 
superior  to  traditional  steady  RANS.  The  methodology  discussed  in  this  paper  has  the  po¬ 
tential  to  make  a  significant  impact  on  these  problems  and  to  bridge  the  gap  between  DNS, 
LES  and  RANS. 
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Figure  1.  Scatterplot  of  the  Yoshizawa  (1986)  model  for  Rmm  versus  the  exact  subgrid  scale 
isotropic  stress  obtained  from  a  96^  direct  numerical  simulation  of  compressible  isotropic 
turbulence  at  a  Mach  number  of  0.1.  (Taken  fcom  Speziale  et  al  1988). 


Figure  2.  Schematic  diagram  of  numerical  simulations  in  turbulence.  (Provided  by  H.  Fasel 
University  of  Arizona). 
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